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ABSTRACT 

We investigate the emission of sub-millimetre-wave radiation from galaxies in a 165 square 
arcminute region surrounding the Hubble Deep Field North. The data were obtained from ded- 
icated observing runs from our group and others using the SCUBA camera on the James Clerk 
Maxwell Telescope, and combined using techniques specifically developed for low signal-to- 
noise source recovery. The resulting 'Super-map' is derived from about 60 shifts of JCMT 
time, taken in a variety of observing modes and chopping strategies, and combined here for 
the first time. At 850 fim we detect 19 sources at > 4ct, including 5 not previously reported. 
We also list an additional 15 sources between 3.5 and 4. Oct (where 2 are expected by chance). 
The 450 /im map contains 5 sources at > 4ct. We present a new estimate of the 850 /im and 
450 fim source counts. The number of sub-mm galaxies we detect account for approximately 
40% of the 850 /im sub-mm background, and we show that mild extrapolations can reproduce 
it entirely. A clustering analysis fails to detect any significant signal in this sample of SCUBA 
detected objects. A companion paper describes the multiwavelength properties of the sources. 



Key words: methods: 
galaxies: formation 



statistical - methods: numerical - large-scale structure of Universe 



1 INTRODUCTION 



The Hubble Deep Field North (HDF-N) <Williams et allll996h . a 
small region of the sky targeted by the Hubble Space Telescope 
(HST), has stimulated the study of the high redshift Universe ever 
since the data were released in 1995. The original optical image of 
the HDF is one of the deepest ever obtained, and resolves thousands 
of galaxies in an area of a few square arcminutes. However, since 
optical images capture only a narrow part of the spectrum, and since 
the rest-frame wavelength range detected depends on the redshift, 
it is necessary to supplement the HST image with data at other 
wavelengths in order to obtain a more complete understanding of 
galaxy evolution. In the years since the HDF image became public, 
deep pointings using radio, X-Ray, near-IR, and mid-IR telescopes 
have been conducted. Additionally, optical spectroscopy has been 
carried out on hundreds of suitable objects in the field, thereby ob- 
taining redshifts for most of the brighter HST-detected objects. 

The original HDF field is the size of a single WFPC2 field 
of view, roughly 2' x 2'. HST also obtained shallower obser- 
vations in fields adjacent to this, extending the region of study 
to roughly 6' x 6'. These 'flanking fields' have also been cov- 
ered by other telescopes, and in fact over time the region associ- 
ated with the HDF has been extended to about 10' x 10' in most 
wave-bands. For an excellent review of the HDF-N region and its 
impac t on the optical view of astronomy, refer to the article by 
iFerguson et alJ l2000l) . Despite these extensive observations, fully 
understanding the high redshift universe is hindered by the pres- 
ence of dust, which re-processes radiation and emits it in the Far- 
Infrared (FIR). The importance of this is clearly demonstrated in 
the population of galaxies detected by the Submillimetre Common 



User Bolometer Array (SCUBA; Holland et al. 1999). Even with 
the tremendous observational effort of the past five years, we have 
made only modest progress with SCUBA detected galaxies beyond 
the co nclus ions drawn from the pioneerin g work of iHughes et alJ 
( 1998) and lSmail. Ivison. & Blair] <1997t) . SCUBA mapping sur- 
veys have constrained the number counts, and show that significant 
evolution is required in the local Ultra-Luminous Infra-Red Galaxy 
(ULIRG) population in order to explain the abundance of high red- 
shift SCUBA sources. 



Having detected these sources, the goal now is to characterise 
them and try to answer fundamental questions about their nature. 
What powers their extreme luminosities? What is their redshift dis- 
tribution? What objects do they correspond to today? One way to do 
this is by doing pointed photometry on a list of high redshift sources 
detected at other wavelengths, but this obviously introduces some 
biases. An alternate approach, one taken by many groups, is to do a 
'blank field' survey, and compare the sub-mm map against images 
taken with other telescopes at other wavelengths. Arguably the best 
field to do this is in the Hubble Deep Field region, which has more 
telescope time invested in observations than any other extragalactic 
deep field. In addition to the data already available, observations 
with the HST-ACS and upcoming confusion-limited SIRTF obser- 
vations, make this an appealing region to target with SCUBA. For 
all these reasons it is worthwhile to combine the available sub-mm 
data in this part of the sky. 
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2 OBSERVATIONS 

Observations with SCUBA require the user to 'chop' the secondary 
mirror at a rate of 7 . 8 1 25 Hz between the target position and a refer- 
ence (also called 'off') position. The difference between the signals 
measured at each position removes common-mode atmospheric 
noise. The direction and size of the chop throw is adjustable, and 
is typically set such that the off position does not change over time 
due to sky rotation. A map made in such a way will then exhibit a 
negative copy of the source in the off position. 

In raster-scan mode, the whole telescope moves on the sky 
in order to sample a region larger than the array size. A map made 
from data collected in this mode is commonly referred to as a 'scan- 
map'. An alternative way to sample a large area is to piece together 
smaller 'jiggle' maps. In this mode, the telescope stares toward the 
target and the secondary mirror steps around a dither pattern de- 
signed to fill in the space on the image plane between the conical 
feedhorns. Photometry mode is very similar to jiggle-mapping, ex- 
cept the dither pattern is smaller in order to spend more time ob- 
serving the target object. The trade off is that the image plane is 
not fully sampled, and there will be gaps in a map made from pho- 
tometry data. Jiggle-mapping and photometry also differ from scan 
mapping because they use two off positions in the chopping pro- 
cess. Therefore there will be two negative echos of each source. 

Our group has produced a shallow 'scan-map' of the region at 
450 and 850 /mi (Borys et al. 2002, hereafter Bo02) using SCUBA. 
Given the high profile of the region however, groups from the 
UK and Hawaii also targeted the HDF to exploit the rich multi- 
wavelength observations available. The ~ 3' x 3' area centred on 
the HDF itself, studied originally by Hughes et al. (1998; hereafter 
H98), was recently re-analysed by Serjeant at al. (2002; hereafter 
S02) and we use their published data for comparison here. Re- 
sults from observations taken by the Hawaii group can be found 
in Barger et al. (2000; hereafter BaOO). We have obtained all of 
these data from the JCMT archive or from the observers directly 
and embedded them in our scan-map. The extra data considerably 
increase the sensitivity of the final map in the overlap regions, but 
at the expense of much more inhomogeneous noise and a corre- 
spondingly more complicated data analysis. Co-adding data taken 
in different observing modes has not previously been performed for 
extra-galactic sub-mm surveys, and so we discuss our procedure in 
some detail below. 

The full list of projects allocated time to study the HDF region 
is extensive. Observing details for each project are given in Table 
Q Almost all work carried out in the region has used the jiggle-map 
mode in deep, yet small surveys. Three projects involved targeted 
photometry observations taken in '2-bolometer chopping' mode. 
The SURF User Reduction Facility (SURF; Jenness & Lightfoot 
1998) software was not able to process these particular files for 
inclusion in the map, but we shall discuss these observations later. 
Project M00BC01 is a single jiggle-map observation taken in such 
a way that co-adding it to the map is also difficult, and shall be 
described later as well. All of the remaining data could be co-added 
into a 'Super-map', as we now describe. 



3 DATA REDUCTION 

We start by simply regridding the bolometer data onto an output 
map, ignoring for a moment the effect of a strongly varying ef- 
fective point-spread function (PSF, i.e. the effect of the different 
chops) across the field. Fig. \\\ demonstrates just how different the 



observations are. Shown is the PSF at various locations in the field, 
alongside the 850 /itn noise map with contours overlaid. Because 
the telescope is moving at a rate of 3 arcsec per sample while tak- 
ing scan-map data, there is no advantage to choosing a pixel size 
smaller than 3 arcsec. Some other studies have chosen smaller (1 
arcsec) pixels for jiggle-maps, but given the uncertainties in JCMT 
pointings, we feel that a choice of 3 arcsecond pixels is sufficient. 

3.1 Flux and pointing calibration 

Registering the individual data-sets with respect to each other is 
difficult, due to the lack of bright sources in the field. The pointing 
log for each night of data was inspected, and in all but a few cases 
there is no reason to distrust the pointing; pointing checks were al- 
ways performed on the same side of the meridian as the HDF, and 
observations were not conducted through transit 1 . Although indi- 
vidual night-to-night pointing comparisons between data sets is not 
feasible, it is still possible to check for a systematic pointing offset 
between scan- and jiggle-maps. One might be concerned about er- 
rors on the order of 3 arcsec, due to the speed at which the telescope 
moves while scanning. 

In a map made solely from the scanning observations, a group 
of sources is detected that correspond to the same objects seen in 
the deep jiggle map in BaOO (project M99AH05). A least squares 
comparison between a 4' x 4' square in the scan-map against the 
same area in the jiggle-map indicates a 4 ± 4 arcsecond offset along 
the direction of the chop (110 degrees east of north). 

Exploiting the well established FIR/radio correlation, we also 
performed a pointing check using the following procedure: We took 
a list of /iJy VLA 1 .4 GHz sources detected by I Richards 2000) and 
extracted the 850 /itn flux from the super-map at each postions. The 
average of these fluxes was used as a figure of merit in checking the 
pointing. By shifting the super-map in both directions and calculat- 
ing this 'stacked' flux, we find that no significant positional offset is 
required maximize the average sub-mm flux at the radio positions. 

For each night of jiggle-map and photometry observations, 
flux calibrations were performed on suitable targets in the area. In 
each case these calibrations were consistent with the averages de- 
termined for others taken in the same time periods, and with the 
larger sample reported in Archibald et al. 1 2002). Since individual 
calibrations may have larger statistical variation, the average cali- 
bration value over the run was used instead. 

Flux calibration of the scan-maps is more problematic. Since 
this mode is currently not so well characterised, it is important to 
compare the map against those taken from the better understood 
jiggle-maps. Although pointing seems well constrained, the fluxes 
in the HDF scan-map are generally larger than their jiggle-map 
counterparts when using the 'standard' gains. This was an issue 
first discussed in Bo02. Unfortunately no usable scan-map calibra- 
tion measurements were taken during the observations. Therefore 
we proceed by adjusting the overall calibration of the scan-map and 
finding, in a chi-squared sense, the value that minimizes the differ- 
ence between the two maps: 

X{ r ) = ^ (^Jiggle - Tuscan) 2 - (!) 

Here the sum is over all the pixels. Because the beam patterns (two- 
versus three-beam) are different, we use only those pixels within a 

1 During much of the period in which this data was obtained, the JCMT 
telescope would drop slightly in elevation while tracking an object through 
transit 
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Table 1. Summary of SCUBA HDF region observations. Type refers to either photometry (P), jiggle-map (J), or scan-map (S). Approximate la noise estimates 
in mjy are provided, though they do vary somewhat across the individual images. The area surveyed in square arcminutes is also listed. Projects denoted by 
an asterisk were photometry observations that could not be added directly to the map, due to the inability of SURF to extract the bolometer positions for 
two-bolometer chopping observations. In addition, M00BC01 was taken using an unfortunate chopping configuration which cancelled out 2 sources near the 
western edge of the map. These data are not folded directly into the map. In total, there have been almost 60 shifts (of 8 hours each) in which the telescope 
was observing this region, which is considerably more than the original HST optical imaging. 



Project ID 


Type 


Shifts 


CT450 (mJy) 


CT850 (mJy) 


Area (arcmin 2 ) 


M97BU65 


P,J 


17 


4 


0.4 


9 


M98AC37 


P 


2 


65 


2 


3 


M98AH06 


P 


4 


11 


1 


3 


M98BC30 


S 


2 


75 


6 


160 


M98BU61 


P 




33 


3 


2 


M99AC29 


S 




95 


6 


160 


M99AH05 


I 


11 


20 


2 


39 


M99BC41* 


P 




N/A 


N/A 


N/A 


M99BC42 


PJ,S 




120 


6 


160 


M00AC16* 


p 




N/A 


N/A 


N/A 


M00AC21* 


p 




N/A 


N/A 


N/A 


M00AH07 


I 


4 


40 


2 


12 


M00BC01 


I 


1 


N/A 


3 


5 


M01AH31 


I 


7 


25 


2 


22 




Figure 1. Layout of the SCUBA HDF observations. The majority of the data is contained within a 12' X 12' square, but other jiggle-maps were also taken 
within the wider Chandra or VLA images. The left panel is a map of the HDF using simulated sources but the real pointing information. This allows the 
variable PSF to be easily seen. Black here corresponds to sources, and white to the negative echoes caused by the necessity to chop. Various parts of the 
field were observed in different ways, as is evident from the different PSFs. On the right is the actual noise map at 850 Lim with contours at 1, 3, 5, 7, and 
9mJy overlaid. The deepest pointings of Hughes et al. 1 1998) are in the centre of the image, while other jiggle-maps and photometry observations can be seen 
superimposed on our scan-map. The corresponding images for the 450 fim map are very similar in appearance, though with higher noise levels. 



beam-width of the sources detected in the jiggle-maps. Based on 
this analysis we determine that scan-map flux conversion factors 
are 0.8 times that of the standard jiggle-map calibrations. 

We obtained scan-map calibration data for other (unrelated to 
HDF) projects that observed planets, and estimated the flux conver- 
sion factor from them. This analysis also obtains the factor of 0.8 
difference between scan- and jiggle-map calibrations. 

3.2 Source detection 

No effort was made to deconvolve this combined map; the observ- 
ing strategies employed by both our group and others do not in- 



terconnect pixels very well, and therefore a robust deconvolution 
cannot be performed. However, in order to gain back the extra sen- 
sitivity from the off-beams, we would like to fit each pixel in the 
map to the multi-beam PSF. This procedure has been adopted by 
other groups as well, but the complication in this case is the vari- 
able PSF across the field. 

Instead of fitting to the PSF, we can fold in the flux from the 
off-beams in a time-wise manner. For each sample, we add the mea- 
sured flux to the pixel being pointed to in addition to its negative 
flux (appropriately weighted) at the position of the off-beam. This 
is equivalent to performing one iteration of the map-making pro- 
cedure described in Johnstone et al. (2000; see also Borys 2002). 
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A single Gaussian is then used to fit for sources in the final map. 
This requires an image relatively free of sources that might lie in 
the location of the off-beam of another source, but produces the 
same output as one would obtain from fitting with the beam pat- 
tern. The resulting image should be considered not so much as a 
map but rather the answer to the question: what is the best estimate 
of the flux of an isolated point source at the position of each pixel? 
The image can only be properly interpreted in conjunction with the 
accompanying noise map. 



4 MONTE-CARLO SIMULATIONS 

There are several simulations that must be performed to assess the 
reliability of the maps. To determine how many detections might 
be false positives, we created a map by replacing the 850 ^im data 
with Gaussian random noise with a variance equivalent to the noise 
estimated for each bolometer. The map was then run through the 
same source-finder algorithm as the real data. This was repeated 
100 times. The average number of positive and negative detections 
as a function of signal-to-noise threshold is plotted in Fig.|2| along 
with the number one would expect based simply on Gaussian statis- 
tics and the number of independent beam-sizes in the map (which is 
an underestimate because it assumes well-behaved noise). A 3.5 a 
cut is commonly used, but these results suggest that at this thresh- 
old about four detections will be spurious in our map (two positive 
and two negative detections). Given that the slope of this plot is 
still rather steep at 3.5 a, small errors in the noise estimate can lead 
to more false positives. Therefore we adopt a 4.0 a cut to deter- 
mine sources in the HDF super-map. Note that the number of false 
positives for the 450 map will be four times larger because the 
beam is half as big (of course if the noise is not well behaved this 
can be even worse!). Therefore, one wants to set a high detection 
threshold for 450 objects (but as we will find later, there is only 
one > 5cr source at 450 /im). The situation is further complicated 
by the fact that the true background is not blank, and is popoulated 
by many unresolved sources that create confusion noise. 

To further estimate the reliability of our detections, we added 
500 sources (one at a time) of known flux for a range of flux levels 
into the map and attempted to extract them using the same pipeline 
as for the real data. A source was considered 'recovered' if it was 
detected with a SNR greater than 4.0 and its position was within 
half the FWHM of the input position. At each input flux, the po- 
sition offset, flux bias, and noise of all the recovered sources were 
averaged and plotted in Figs.[3]and|4]for the 450 /im and 850 pim 
maps, respectively. 

The panel showing completeness is self explanatory; it plots 
the percentage of sources recovered in the Monte-Carlos as a func- 
tion of input flux. As one expects, it is difficult to recover faint 
sources around the noise limit, but very bright sources are always 
recovered. We will discuss this issue more when deriving source- 
counts. 

The ratio of output and input flux in the adjacent panel shows 
that sources fainter than the threshold have been scattered up due 
to the presence of noise and are 'detected' with higher than their 
true flux. The relevant noise components in this case are instru- 
ment al noise and con fusion noise. The latter is called Eddington 
bias lEddingtonll 1 94fJ) . and dominates the central part of the map 
where the instrumental noise is smallest. At 850 fim, the confusion 
limit is ~ 1 mjy, and thus sources with instrumental noises around 
this level will be subject to Eddington bias. 

The RMS of the difference between input and output positions 



i 1 1 1 1 i 1 1 1 1 i 1 1 1 1 r 




3.0 3.5 4.0 4.5 

a threshold 



Figure 2. Number of positive sources in the 850 (im map expected at ran- 
dom. The filled circles show how many spurious sources were detected in 
simulated maps made from noise alone. Error bars are omitted for clarity, 
but are small because of the large number of sources generated in the many 
simulations. The number of positive sources one would expect by chance 
given Gaussian statistics and a Gaussian-shaped beam with a FWHM of 
14". 7 in an 165 square arcminute region is plotted as the thick line. 



allows us to estimate a positional error for real detections as a func- 
tion of flux level. As one would expect, the uncertainty is smaller 
as the input flux goes up. The final plot in the sequence shows the 
average noise level associated with the recovered sources as a func- 
tion of input flux. At the faint end, the noise is lower because such 
sources are only detected in the deepest parts of the map. As the 
flux of the source increases above the noise level of the least sen- 
sitive region (in this case the underlying scan-map), the average 
noise of the detections levels off to the average noise level of the 
field. The 450 /j,m plots (Fig. 3) are very similar to those at 850 /j,m 
(Fig. 4), scaled by approximately a factor of 10 in flux. 

Basic conclusions at 850 fim are that at a flux limit of 8 mjy 
the source counts are about 80% complete, fluxes are biased by 
only a few percent, and positions are accurate to about 3.5 arcsec. 
The brighter objects are constrained much better than this, but there 
are few of them. At the faintest levels confusion has a significant 
affect on fluxes, positions, and completeness. 



5 SUB-MM SOURCES IN THE HDF 

Most of the 850 /im sources exhibit off-beam signatures that are 
distinguishable by eye (the negative echos to the left and right of 
the sources). Finding sources at 450 /im, however is more difficult. 
The single beam pattern is not well described by a Gaussian, plus 
the sensitivity at 450 ^im is too poor to detect any but the bright- 
est sources. In addition, being more weather-dependent, the noise 
is even more inhomogeneous than for the 850 /im data. Neverthe- 
less, we will report the 3<j upper limit to the 450 /im flux for each 
850 fim detection. To avoid reporting spurious detections, we have 
set a threshold of 4.0cr on the 850 /im catalogue derived from the 
super-map. However, a supplementary list of sources at 850 /j,m de- 
tected above 3.5ct is also provided for comparison against other 
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Table 2. Sub-mm detections in the HDF super-map. Only sources with a SNR (shown in brackets) greater than 3.5 are listed. If the source was previous! 
detect ed in a separate sur vey, t hen the reference to the paper is given, as well as the measured flux from there. The papers are BaOO: lBarger et aT . 2000; 
Bo02: <Borvs et all20 02bl S02: lSerieant et all2002l) . 3cr upper limits in the other channel are given for each source. 



ID 


RA 


DEC 


5850 ± °"850 


5450 ± f450 


Previously Detected 


850 fj,m detections > 4(7 


SMMJ123607+621145 


12:36:07.3 


62:11:45 


15.2 ±3.8 (4.0ct) 


< 97 


Bo02: 15.4 ± 3.4 mJy 


SMMJ123608+621251 


12:36:08.5 


62:12:51 


16.0 ±3.7(4.3<t) 


< 95 


Bo02: 13.8±3.3mJy 


SMMJ 123616+6215 18 


12:36:16.6 


62:15:18 


6.3 ±0.9 (6.5(7) 


< 35 




SMMJ123618+621009 


12:36:18.9 


62:10:09 


6.6 ± 1.5(4.3(7) 


< 66 




SMMJ123618+621554 


12:36:18.7 


62:15:54 


7.2 ±0.9 (8.1(7) 


< 34 


BaOO: 7.8±1.6mjy 


SMMJ123621+621254 


12:36:21.8 


62:12:54 


12.1 ±2.6(4.7cr) 


< 85 


Bo02: 11.4 ± 2.8mJy 


SMMJ123621+621712 


12:36:21.3 


62:17:12 


8.8 ± 1.5(6.0(7) 


< 72 


BaOO: 7.5 ± 2.3 mJy 












Bo02: 13.2 ± 2.9 mJy 


SMMJ 123622+621618 


12:36:22.6 


62:16:18 


8.6 ± 1.0(8.5(7) 


< 45 


BaOO: 7.1 ± 1.7mJy 


SMMJ123634+621409 


12:36:34.2 


62:14:09 


11.2 ± 1.6(7.1cr) 


< 67 


BaOO: 19.0±2.8mJy 


SMMJ 1 23637+62 1157 


12:36:37.3 


62:11:57 


7.0 ±0.8 (8.3(7) 


< 46 




SMMJ123646+621451 


12:36:46.3 


62:14:51 


8.5 ± 1.3(6.4(7) 


< 47 


BaOO: 10.7 ± 2.1 mJy 












Bo02: 11.4 ± 2.9mJy 


SMMJ123650+621318 


12:36:50.6 


62:13:18 


2.0 ±0.4 (5.3(7) 


< 11 


S02: HDF850.4/5 2.1 ± 0.3 mJy 


SMMJ123652+621227 


12:36:52.3 


62:12:27 


5.9 ±0.3 (18.0(7) 


< 12 


S02: HDF850.1 5.6 ± 0.4mJy 


SMMJ 123656+62 1203 


12:36:56.6 


62:12:03 


3.7 ±0.4(8.8(7) 


< 16 


S02: HDF850.2 3.5 ± 0.5 mJy 


SMMJ 123700+6209 12 


12:37:00.4 


62:09:12 


8.6 ±2.1 (4.1(7) 


< 85 


BaOO: 11.9±3.0mJy 


SMMJ123701+621148 


12:37:01.3 


62:11:48 


4.0 ±0.8 (5.2ct) 


< 27 


S02: HDF850.6 6.4 ± 1.7mJy 


SMMJ 123703+62 1303 


12:37:03.0 


62:13:03 


3.4 ±0.6 (5.3ct) 


< 25 




SMMJ123707+621412 


12:37:07.7 


62:14:12 


9.9 ± 2.5 (4.0(7) 


< 85 




SMMJ123713+621206 


12:37:13.3 


62:12:06 


6.1 ± 1.4(4.3(7) 


< 40 


BaOO: 8.8±2.0mjy 


Additional 850 fim detections > 3.5(7 and < 4. Oct 


SMMJ123607+621021 


12:36:07.3 


62:10:21 


13.5 ±3.7(3.7cr) 


< 97 




SMMJ123608+621433 


12:36:08.5 


62:14:33 


6.1 ± 1.7(3.6(7) 


< 60 




SMMJ12361 1+621215 


12:36:11.9 


62:12:15 


12.8 ±3.4 (3.7cr) 


< 90 


Bo02: 12.2±3.0mJy 


SMMJ123628+621048 


12:36:28.7 


62:10:48 


4.4 ± 1.2(3.7(7) 


< 55 




SMMJ123635+621239 


12:36:35.6 


62:12:39 


3.0 ±0.8(3.6a) 


< 41 


S02: HDF850.7 5.5 ± 1.5 mJy 


SMMJ 123636+620700 


12:36:36.9 


62:07:00 


22.1 ±5.6 (3.9(7) 


< 155 




SMMJ 123648+62 1842 


12:36:48.0 


62:18:42 


19.5 ±5.4 (3.6(7) 


< 154 




SMMJ 123652+62 1354 


12:36:52.7 


62:13:54 


1.8 ±0.4(3.9a) 


< 16 


S02: HDF850.8 1.7 ± 0.5 mJy 


SMMJ123653+621121 


12:36:53.1 


62:11:21 


2.8 ±0.8(3.6(7) 


< 40 




SMMJ123659+621454 


12:36:59.1 


62:14:54 


5.2 ± 1.4(3.8(7) 


< 72 




SMMJ 123706+62 1851 


12:37:06.9 


62:18:51 


21.6 ±5.8 (3.8(7) 


< 178 




SMMJ123719+621109 


12:37:19.7 


62:11:09 


7.2 ± 2.0(3.6(7) 


< 55 




SMMJ123730+621057 


12:37:30.4 


62:10:57 


13.3 ±3.6 (3.7(7) 


< 98 


Bo02: 14.3 ± 3.2 mJy 


SMMJ 12373 1+62 1857 


12:37:31.0 


62:18:57 


27.1 ±7.6(3.6(7) 


< 286 




SMMJ123741+621227 


12:37:41.6 


62:12:27 


23.7 ±6.1 (3.9(7) 


< 185 




450 (im detections > 4cr 


SMMJ 1236 19+62 1127 


12:36:19.3 


62:11:27 


< 5.8 


110 ± 26 (4.2cr) 




SMMJ 123632+62 1542 


12:36:32.9 


62:15:42 


< 5.9 


105 ± 25 (4.2cr) 




SMMJ 123702+62 1009 


12:37:02.6 


62:10:09 


< 5.2 


120 ± 27(4.4cr) 




SMMJ123727+621042 


12:37:27.4 


62:10:42 


< 10.4 


220 ± 42 (5.2cr) 




SMMJ 123743+62 1609 


12:37:43.8 


62:16:09 


< 24.0 


300 ± 72 (4.2cr) 





data sets. The full list of 850 fim and 450 fim sources is presented 
in Table|2| 

A simple test of source reality is to search for negative flux ob- 
jects. All but two negative sources found were associated with the 
off-beam of one of our positive detections. Based on the Monte- 
Carlos, this is not unreasonable. There are 19 sources at 850 
detected over 4cr, 5 of which are new, not having been reported be- 
fore in the individual surveys. Apart from a few exceptional cases, 
described below, all sources previously reported in the region are 
recovered at comparable flux levels. An additional 15 sources are 
detected between 3.5a and 4cr. Our Monte-Carlos suggest that on 



average 2 of these may be spurious. An image of the 850 fim super- 
map is given in Fig. [5] There are only 5 sources recovered from 
the 450 fim image, which is shown in Fig. [6] There have been no 
450 fim detections previously reported in the HDF. The finder chart 
in Fig.Qcan be used to identify the sources extracted by our algo- 
rithm in each of the maps. 

5.1 Comparing the source list against previous surveys 

A plot comparing the recovered fluxes against those previously 
published is given in Fig. [8] We show the difference between the 
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Figure 5. The 850 fj,m signal-to-noise map. The greyscale is stretched to highlight sources. This image can be used in conjunction with the finder chart in 
Fig.Qto identify the objects reported in the source list. 



published fluxes and super-map estimates. The flux differences 
should be zero, and based on the size of the error bars it is clear 
that no significant variations exist between the estimates. The ob- 
jects that appear slightly discrepant are discussed at length below, 
and all sources will be described in detail in a forthcoming paper. 
It should be noted that the error bars on the sources recovered from 
the super-map are generally smaller than those obtained from the 
individual sub-maps. 

5.1.1 The central HDF region from Hughes/Serjeant et al. 

Eight sources, labeled HDF850. 1 through HDF850.8 are associated 
with the data collected by H98. The more thorough analysis in S02 
found that one of them (HDF850.3) detected in the original map 
fails to meet new detection criteria. It is also the only source unde- 
tected in the super-map presented here. HDF850.4 and HDF850.5 
appeared to be a blend of two sources in their original map, and 
therefore both H98 and S02 took the extra step of attempting to fit 



the amplitude and position of them. With the super-map however, 
this pair of sources is better fit by a single source with flux compa- 
rable to the sum of fluxes from the two extracted by S02. This may 
be because we use larger pixels (3" compared to 1" from S02) and 
therefore lose some of the resolution required to separate adjacent 
objects. Sources HDF850.6 and HDF850.7 are detected at a fainter 
level than reported in S02. These discrepancies might be partly due 
to their positions - both are in a noisier region of the individual sub- 
map, and near its edge. The super-map contains additional data (es- 
pecially photometry) in the central HDF region, and so one would 
expect our noise estimates to be mildly lower compared to those of 
S02, which is indeed generally true. 

5.1.2 HDF flanking field jiggle-maps from Barger et al. 

The seven sources detected individually in BaOO are confirmed in 
the super-map at comparable flux levels except for one object not 
detected in the radio. This is most likely due to differences in anal- 
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ysis methods; BaOO used aperture photometry with annuli centred 
on radio positions believed to be associated with the sub-mm detec- 
tion. The source in question had no counterpart, and therefore de- 
termining flux with aperture photometry is not as straightforward. 
Even when the sub-map alone is considered individually, the mea- 
sured flux does not match that originally reported. 



5.1.3 Scan-map observations ofBorys et al. 

All six sources from the 4a list of Bo02 are recovered, although two 
of them are detected at lower significance in the super-map. Four 
of the six sources from the supplementary list of 3.5a sources are 
not recovered in the super-map. Three of these exist in regions of 
overlap between surveys. This supports the warning made in Bo02 
that 4a SCUBA detections are less likely to be spurious than lower 
SNR ones. In general, however, we find that sources detected at 
> 4a are confirmed in separate sub-maps. This is an important test 
of the reliability of faint SCUBA detections, which can only be ef- 



fectively carried out in the HDF region, because of the overlapping 
independent data-sets. 

5.2 Comparison with photometry observations 

As we have mentioned, a number of photometry observations taken 
in the two-bolometer chopping mode were conducted. Although 
these data sets cannot be co-added into the map, we can compare 
flux estimates for the photometry bolometer with the position in the 
super-map. 

We should also note that while the scan-map observations 
were being taken, several photometry observations were conducted 
at positions of tentative detections in order to verify them. None 
of these positions correspond to detections in the final map, but 
all have fluxes consistent with the corresponding position in the 
super-map. This illustrates that the practice of picking out low SNR 
sources 'by eye' in SCUBA surveys is not very effective. Our group 
has also conducted two observing programmes designed to un- 
derstand the sub-mm properties of Lyman break galaxies (LBGs, 
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Figure 7. Finder chart for sub-mm sources detected in the HDF super-map. The outlines a re 850 Aim noise conto ur levels from 1 to 1 1 mJy (darkest) in steps 
of 2mjy. Tri angles are 850 lira s ource positions from lBarger et alJ 1*2000)1. diamonds from Serj eant et all squares with an embedded plus are the 4<r 

sources from lBorvs et alj<2002bl) . and the plus symbols are the lower SNR sources from that survey. Circles with crosses are sources recovered above 4<r from 
the super-map, and crosses show those additional sources detected at lower significance. The five 450 sources detected in the super-map are denoted with 
asterisks. 



Chapman et al. 2000) and optically faint radio sources (OFRS, 
Chapman et al. 2002). These types of source will be discussed fur- 
ther in paper II. Three measured LBGs that fall within the region 
of the HDF are not detected in either the photometry measurements 
or the super-map. Of the seven OFRS photometry measurements, 
all have comparable fluxes at their corresponding position in the 
super-map. A summary plot similar to the one shown previously 
when comparing known detections is given in Fig.|9| Again, no 
significant discrepancies exist between the estimates. 



6 NUMBER COUNTS OF SUB-MM SOURCES 

In order to estimate the density of sources brighter than some flux 
threshold S, N(>S) from our list of detected sources in Table|2| 
we must account for several anticipated statistical effects: 

(i) The threshold for source detection, St = rna is not uniform 
because the noise varies dramatically across the map. 

(ii) Due to confusion and detector noise, sources dimmer than 
St might be scattered above the detection threshold and claimed 
as detected. 

(iii) Similarly, sources brighter than St might be missed be- 
cause of edge effects, possible source overlaps, and confusion. 



Item (iii) is simply the completeness of our list of sources. 
For a source density N(S)dS which falls with increasing flux, the 
effect of item (ii) can exceed that from item (iii), resulting in an 
Eddington bias in the estimated source counts. One can calculate 
the ratio of the integrated source count to the number of sources 
detected using the 'detectability': 

1 J™N(S)dS 
j(S') J-^S^NWdS' 

where N(S)dS is the number of sources with a flux between S 
and S + dS. We have introduced the quantity <j)(S, S') which is the 
fraction of sources between a flux, S and S + dS, that are detected 
above a threshold, S'. Note that (j)(S, S') ranges between zero and 
unity, but 1/7(S') can be larger than one, depending on the form of 
N(S) and choice of S'. For bright sources, where the completeness 
is 1.0 and the SNR is high, 0(S, S') should approach unity. 

The numerator of equation|2|is the quantity we are trying to 
determine, and the denominator is the output from the survey. To 
determine N(> S) we simply take the raw counts from our survey 
and multiply them by 7(S'). The calculation of 7(S') from the 
Monte-Carlo estimates of <f>(S, S') requires a model of the source 
counts, of which there are many different forms in the literature. All 
share the property that they are steep at the bright end and shallow 
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Figure 3. 450 /im source recovery Monte-Carlo results. From left to right 
and top to bottom the plots are: a) completeness; b) flux bias ratio; c) point- 
ing error; and d) noise level. Error bars are not plotted for clarity, but one 
can gauge the uncertainty by the scatter of the points in the vertical direc- 
tion. 



Serjeant Barger 



Borys 




Figure 8. Our new source list compared with previously reported detec- 
tions. We show the difference between previously reported flux and that 
derived here for each source in common. The outline denotes the size of the 
error-bars from the previously reported estimates. We subtract the reported 
fluxes from our own, and plot the difference and error-bar as the filled trian- 
gles. See the text for a more complete description for those few points that 
seem discordant. 
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Figure 4. 850 /im source recovery Monte-Carlo results. The order is the 
same as in the last figure. In general the plots are very similar to their 
450 (im counterparts, with fluxes an order of magnitude smaller. 




Figure 9. Comparing super-map fluxes with photometry estimates. We plot 
the difference between reported photometry fluxes from observations not 
added into the map and those derived from the super-map. As before, the 
outline shows the error-bar from the reported values, while our points (after 
subtracting off the reported fluxes) are plotted as solid triangles. Since some 
of the photometry goes much deeper than our map, their error estimates are 
in those cases much lower. All sources are in reasonable agreement within 
the combined error-bar of each estimate. Here, 'LBG' is a Lyman-Break 
Galaxy target and 'OFRS' is an optically faint radio source. 
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at the faint end. They are also typically constrained such that the 
total amount of energy does not exceed the measured value of the 
sub-mm extra-galactic background. We employ the double-power 
law form described in IScottetaljEool . 



dN(> S) 
dS 



No 
So 



(!)"+( 



So) 



(3) 



and use So = 1.8 mjy, a = 1.0, (3 = 3.3, and 7V = 1.5 x 
10 4 deg" 2 . 

Our Monte-Carlos give us <f)(S,S'), which we fit using the 

form 

<p(S, S') = l- exp (A(S - B) c ) . (4) 

for each value of S". This was chosen empirically to match the 
shape of <f)(S,S'). An alternative approach would be to perform 
more Monte-Carlos with a finer spacing in AS and then spline the 
result to allow for interpolation between the sample points - the 
results are almost identical. Using these estimates of <f>(S, S'), to- 
gether with the source-count model, ^y(S') can be computed nu- 
merically using equation|2| 

The quantity 7 relates what we want to determine, N(>S), 
with the number of sources our survey detects. Obviously 7 is in- 
fluenced by what model is used in the calculation described above, 
and other forms of the source spectrum are found in the literature. 
We find 7 varies by no more than 15% across a wide range of rea- 
sonable parameter values. This is smaller than the Poissonian error 
caused by having so few sources detected. 



6.1 The 850 /j,m number counts 

In Fig. 1101 we plot the 7 determined using the source count pa- 
rameters described above, as well as <j}(S, S') determined from the 
Monte-Carlos. 

The detectability passes through unity at around 7mJy; 
sources fainter than this are not detected very efficiently, and there- 
fore we must boost the raw count to account for the incomplete- 
ness. Past this, brighter source counts get a slight boost from fainter 
sources that have been scattered above the threshold due to noise. 
Given 7(S), we can now calculate the counts based on the number 
of sources detected in the entire 165 square arcminute super-map. 
This is given in tabular form (see Table|5J and as a plot alongside 
other estimates of the 850 /im source counts (in Fig. l 1 1 i . 

For the above source count approach we have provided esti- 
mates at a number of flux values (chosen to be integer numbers of 
mjy). There are other ways of presenting the counts which avoid 
the need for binning. An alternative method uses the flux bias ra- 
tio (e.g. a smooth curve fitted to the bottom right panel of Fig. 4) to 
transform the flux of each detected source, and then estimate the ef- 
fective area in which such a source could have been detected. The 
results of this approach are shown by the step function plotted in 
Fig. ll II and are clearly consistent with the other estimate. 

In general the derived counts are in excellent agreement with 
estimates from other surveys. They also agree broadly with the 
input model used to estimate 7(5). It should be noted that while 
other surveys use 68% confidence bounds, we prefer to quote 95% 
limits. It is interesting to note that the bright counts from this sur- 
vey and the UK 8-mJy project are slightly lower than those found 
from surveys of smal ler areas w hich often involve the mo deling 
of cluster lenses jCowie. Barger. & Kneibl2 002: Smai l et all2002t 
Chap man et all2 0021. Comparing the ratio of counts derived from 
our survey relative to the 8-mJy project yields a factor of ~ 0.85 




S/mJy 

Figure 10. Summary of the 850 fj,m source-count calculation. In the top 
panel we plot the detectability, determined using the source-count model 
described in the text. As expected, it is higher at low flux levels (indicating 
that the survey has missed sources) and approaches unity at the bright end. 
Between 7 and 12mJy the detectability is slightly lower than one, demon- 
strating that some sources below the threshold are contributing. The error 
bars represent the scatter obtained when using a range of source-count mod- 
els that roughly fit the measured counts and sub-mm background. In the 
lower panel, we plot ij>(S,S') for 5" ranging from 0-16 mjy in steps of 
4 mjy. The S' = case is simply the completeness of the survey. 



for the 6-12 mjy range. If we fit to the average of the cluster counts 
and interpolate across this same range, this ratio between our counts 
and theirs is roughly 0.4. This might indicate that those surveys had 
been sensitive to clustering, which would naturally lead to an over- 
estimate of the number counts. It could also be that the cluster lens 
surveys were contaminated by sources intrinsic to the clusters, or 
that there was some systematic bias in the lens models. We also 
note that some of the surveys include candidates detected with a 
significance lower than the 4<r cut we used here. 

Although the current survey is limited by the lack of bright 
sources in the HDF, it is clear that the counts fall off quite steeply 
with increasing flux. Different surveys are often compared via the 
slope of a power-law form of the number counts (N(S)dS oc 
S~ a ). For our counts we obtain a — 2.9 ± 0.5 (95% confidence 
limits), which is in excellent agreem ent wi t h estim ates obtain by 
other groups: 2.8 ± 0.7 from iBlain et alj <1999l) . 3. 2+°I from 
IB areer. Cowie. and S andersl 1 1 999h . and 3 . 2 ± 0.7 from Eales et al. 
feOOOh . 

At the faint end, it seems that the counts turn over and flat- 
ten out with decreasing flux. Indeed if they do not then the sub- 
mm background will be overproduced. Deeper surveys that reach 
to lower flux levels will be needed in order to determine this unam- 
biguously. 



6.2 The 450 ^tm number counts 

Since we have only detected 5 objects at 450 [im, we choose to 
quote only a single value for N(> 100 mjy). Because the num- 
ber counts are not well constrained at all at 450 /J,m it is difficult 
to know what model to use in the completeness simulations. We 
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Table 3. 850 fim source counts from the HDF super-map using only the > 4cr detections. The total area surveyed is 165 arcmin 2 , and the error ranges are 
95% confidence bounds on a Poisson distribution. We list results in steps of 2 mjy. We also include published results from other surveys for comparison. These 
include the comparable results from the wide area 8-mJy blank field survey, and three separate surveys conducted toward galaxy clusters, which appear to 
show a higher estimate than blank field counts. 



S[mJy] iV raw (> S) 7(5) JV(> S) 



JVgmJyO S) 



^cluster (> S) 

JCowi^aige^«nei^00Z| 



^cluster (> S) 
(Smai] ci al. 2002) 
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Table 4. 450 (im source counts from the HDF super-map. We present here 
our estimate and those compiled from the literature. The other estimates are 
based on a different approach and also chose to quote symmetric error bars. 



Flux (mjy) 


N(> S) 


Comment 


10 

25 
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2100+ 1200 
ziuu_ 1200 

OUU -500 
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Smail et al. (20021 
Smail et al. (2002) 
This work 



chose a simple power law form with N(S)dS oc S~ a . Our esti- 
mates of the detectability were largely insensitive to the choice of 
a, and turned out to be simply the inverse of the completeness at 
100 mjy (80%). This result, along with two previous estimates at 
lower fluxes, is given in Table|4| There are no published estimates 
at the bright flux levels which are probed here. 

A fit to a power law with these combined counts gives a = 
2.2 ± 0.4. This is not too dissimilar from the 850 fim slope, 
though slightly shallower. This might be an interesting result: if 
the SCUBA sources are mainly at z ^ 1, then one would expect 
the 450 fim fluxes to drop off more steeply than the 850 counts 
because the K-correction becomes positive for 450 /im at these 
redshifts. With only three estimates for the 450 /im number counts 
however, it is premature to draw strong conclusions. Nevertheless, 
this shows that constraints on the number counts at different wave- 
lengths can serve as a probe of the evolution and redshift distribu- 
tion of sub-mm galaxies. 



7 CLUSTERING OF SUB-MM SOURCES 

As we have already mentioned, the pr esence of clustering m ay also 
affect estimates for the source counts. IPeacock et alj i2000) report 
weak evidence of clustering in the ~ 2 x 2 arcminute map of H98 
in the sense of statistical correlations with LBGs, which are them- 
selves strongl y clustered jGiavali sco & Dickinsonll200ll) . The UK 
8mJy survey IScott et alj |2002) which covers over 250 square ar- 



cminutes, and the smaller yet deeper CFRS 3- and 14-hour fields 
(Web b et all2003l) fail to detect clustering among the SCUBA de- 
tected sources. There appears (by eye) to be clustering in the super- 
map; in particular the concentration of sources near the centre of 
the map, the trio of sources to the west of the map, and the group of 
4 north of it might suggest a clustering scale on the order of 30 arc- 
sec or so. However, one would expect more sources in these areas 
because of the increased sensitivity. Hence one needs to carry out a 
statistical clustering analysis, including the inhomogeneous noise, 
to quantify this. 

Clustering is usually described as the probability, p, of finding 
a source in a solid angle fl, and another object in another solid 
angle Q separated by an angle 8. This probability is described by: 



p(0) = N 2 [l + w(6)}Q 2 



(5) 



where N is the mean surface density of objects on the sky and w(6) 
is the angular two-point correlation function. If w(8) is zero, then 
the distribution of sources is completely random, while otherwise 
it describes the probability in excess of random. The correlation 
function w(0) is estimated by counting pairs and there are several 
specific estimators in the literature. The one we employ is that pro- 
posed bv lLandv and Szalavlil99l : 



w{8) 
with Sw(8) 



DD - 2DR + RR 
RR ' 



■w(6) 



DD 



(6) 
(7) 



This particular estimator has been shown to have no bias and also 
has a lower variance than the alternatives. In this equation, D repre- 
sents sources in the SCUBA catalogue, and R are sources recovered 
from Monte-Carlo catalogues. DD is the number of pairs of real 
sources that fall within a bin of width 88 in the map. DR are data- 
random pairs, and RR are random-random. For simplicity, each 
catalogue is normalised to have the same number of objects. To ob- 
tain the random catalogues, we created 1000 mock fields based on 
the source count model used in the previous section, and the actual 
noise of the real data. The sources were placed randomly through- 
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Figure 11. The 850 (im cumulative source counts. The solid diamonds denote the results from the c urrent wo rk. The step-function sh ows counts de rived from 
the same data, but using a different recipe (see Section 6.1). Counts derived from cluster studi es by Chapman et al. 1 2002) and Smail et al. 1 2002) are shown 
by the open squares and circles respectively. The U K 8-mJy survey cou nts IScott et aT] 120021) are shown as open triangles. Stars represent the counts from 
iHughes et alj Il998t) - and the open diamond is from lBorvs et all [2002b). Some points are slightly offset along the flux axis for clarity. Overlaid is the two 
power-law model used in the calculation of 7 (solid line), and two predictions based on representative galaxy evolution models from Rowan-Robinson 1200 ll) 
- the dashed line represents a universe with Qyi = 1.0 and Ha = 0.0 while the dotted line is Qm = 0.3 and Ha = 0.7. The dot-dashed line below the 
curves is the count prediction obtained from extrapolating the IRAS 60 /im counts and invoking no evolution. 



out the field, and the resulting mock data were placed into the same 
pipeline as our real data. 

This a pproach is s lightly different from that of I Webb et alJ 

(2003) and IScott et alJ 120021) . the only other two surveys to at- 
tempt a clustering analysis of SCUBA sources. In those analyses, 
the mock images were modified only by adding noise to each pixel. 
The amount of noise added was taken from the noise map which 
was created along with the real signal map. Therefore their final 
mock images do not exhibit the chop pa ttern that one wou ld ex- 
pect to see. Recognizing this limitation, Web b et alJ I <2003l) took 
the added step of masking out regions in the mock images that 
correspond to the positions of the off-beams in the real map. Our 
simulations involve a full sampling of the mock images using the 



astrometry information from the real data. Therefore the simulated 
and real maps have the same beam features. 

We used relatively wide 30 arcsec bins in because the num- 
ber of sources is so low. This is also twice the size of the SCUBA 
beam at 850 fim, which was the binsize criterion used bv lScott et alJ 
(2002). w(6) is estimated using both the 3.5 and 4<r catalogues. 
Though some of 3.5a sources may be spurious, the increased num- 
ber of objects helps bring down the clustering error bars. As Fig. 1121 
shows, however, even when the 3.5a sources are included there is 
no evidence for clustering in the HDF super-map, since there is no 
angular bin that has a w(6) significantly different from zero. For 
comparison in that figur e, we also plot esti mates of w(d) for EROs 
and LBGs foaddi et alj|200fi iGiavalisco & Dickinsonll200ll) . We 
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Figure 12. Angular two-point correlation function estimate for the HDF 
super-map sources. Our estimates of w(8) are shown as the open circles 
(for the > 3.5<t catalogue) and filled circles (for the > 4. Oct catalogue). 
The points are plotted at the midpoint of each 30 arcsec bin. The first point 
is zero because there are no objects in the super-map closer than 30 arcsec 
to another. The measured clustering signal of EROs is shown as a dashed 
line, and that of LBGs as a dotted line. There is no evidence of sub-mm 
clustering here, though the errors are still quite large. 



repeated this using bins half as wide, and again no significant devi- 
ation from zero was found. 

This is not the only clustering estimator one can use. We also 
performed a 'nearest-neighbour' analysis (see e.g. Scott & Tout 
1989) to test if sources were closer together than expected at ran- 
dom. This statistic simply examines the distribution of distances to 
the nearest neighbours for each source. The cummulant of nearest- 
neighbours is then compared against our set of Monte-Carlo cata- 
logues to determine if a pair-wise clustering signal is present. The 
results are shown in Fig. 1131 

There is a lack of sources with neighbors closer than ~ 30 
arcseconds, but a formal Kolmogorov-Smirnov test indicates only 
a ~ 40% chance that the distributions are significantly different. 
So there is no strong evidence of clustering here either. 

Note that in each of these clustering analyses, it is difficult to 
estimate the clustering strength on scales near to the beam size. Our 
source extraction algorithm is insensitive to a fainter source closer 
than 12 arcsec to a brighter one. Note that if SCUBA sources are 
clustered in a similar manner as EROs or LBGs, then the signal 
in w(6) would be expected to be strongest at <C, 30 arcsec. This 
is only a factor of 2 larger than the SCUBA beam, and therefore 
one expects that a clustering detection with SCUBA will be rather 
difficult due to blending of the sources. Our particular approach 
has not been optimised for separating nearby sources and hence 
our catalogue is not ideal for clustering analysis. Future attempts to 
measure the clustering strength should pay particular care to source 
extraction algorithms. 

There is currently no convincing detection of sub-mm clus- 
tering; what is needed is a very large (~ 1 square degree) sur- 
vey with hundreds of sources in order to decrease the error bars 
on the clustering estimate. The on-going SHADES survey (see 
http://www.roe.ac.uk/ifa/shades/ ) purports to do just that. 



Figure 13. Nearest neighbour clustering analysis. In the top panel we plot 
the nearest neighbour distribution of the data compared with a Monte-Carlo 
generated random set. The random data falls off below 12 arcseconds, 
where our source extractor is unable to seperate unique objects. There is 
a lack of real sources detected with pair separations less than ~ 30 arcsec. 
This is also reflected in the bottom panel, where we plot the cumulative dis- 
tribution of nearest neighbours (dots) compared with Monte-Carlos (solid 
line). 



8 IMPLICATIONS OF THE COUNTS 

8.1 The 850 /im sub-mm background 

We explored the range of parameters for the two-parameter model 
that could fit our source counts and still fall within the limits of the 
FIR 850 jim background constraint: 3.1 — 4.1 x 10 4 mj y deg -2 
IPuget et all 199d Mauser et alJll998HFixsen et alll998h . By cal- 
culating the integral of SN(S)dS we find that our > 2 mJy sources 
contribute 1.4 x 10 4 mjy deg -2 to the FIB. This is consistent with 
estimates from several groups, and demonstrates that a significant 
fraction of the sub-mm Universe is still below the flux threshold at- 
tainable from current SCUBA surveys. However, given the freedom 
which still exists in the faint end counts, and in the current level of 
uncertainty of the sub-mm background itself, the entire background 
can easily be made up of sources with Ssso^m > 0.1 mjy. This re- 
sult can be used to constrain models that predict the evolution of 
IR luminous galaxies. Future surveys that detect more sources will 
constrain the source counts further, and extend the limiting flux 
down to fainter levels. 

8.2 Evolution of sub-mm sources 

The counts are much greater than what one obtains by estimating 
85 fjm fluxes simply from the IRAS 60 /im counts. As addressed 
bv lBlain et at] J2002I) . a SCUBA galaxy with a flux ^ 5 mjy has 
an inferred luminosity in excess of 10 12 L Q if they are distributed 
at redshifts greater than 0.5. Note that Chapman et al. (2003) find 
spectroscopic redshifts z > 0.8 for each source in their sample of 
SCUBA detected radio galaxies. Such sources make up at least 50% 
of the SCUBA population brighter than 5 mjy, and radio undetected 
sources are likely to be at even higher redshifts. 

At these luminosities the number of SCUBA objects per co- 
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moving volume is several hundred times greater than it is today. 
Therefore there must be significant evolution past z J> 0.3 (the 
IRAS limiting redshift). 

Modeling this evolution has been difficult due to the lack of 
observational data on the redshift distribution of SCUBA sources. 
Attempts to model the luminosity evolution of the SCUBA sources 
have been carried out using semi-analytic methods (e.g. Guiderdoni 
et al. 1998) and parametric ones (e.g. Blain et al. 1999, Rowan- 
Robinson 2001, Chapman et al. 2003). In many cases, the starting 
point is the well determined IRAS luminosity function. This gives 
us the number of sources of a given luminosity per co-moving vol- 
ume. These luminosities are then modified as a function of redshift, 
and then 850 /im fluxes are extrapolated and source counts deter- 
mined. We expect the number of galaxies to increase in the past 
due to merger activity, so number evolution must play some role, 
but it is noted that strong number evolution overproduces the FIR 
background. There are many models (and a range of parameters 
within them) that fit the current data. Therefore, until we can better 
constrain the counts and determine redshifts, the only firm conclu- 
sion one can make is that SCUBA sources do evolve strongly. Our 
catalogue of 34 SCUBA sources within the HDF region, with the 
amount of multi-wavelength data being collected there, should con- 
tribute toward distinguishing between models. 



8.3 A connection with modern-day elliptical galaxies? 

Observations of local star-forming galaxies suggest that 

SFR = 2.1 x 10" 10 (L 6OMm /L ) M yr-\ (8) 

where Leo/im = (z^SV) leo^m is me 60 ^im luminosity (Rowan- 
Robinson et al. 1997). Assuming the SLUGS result for dust tem- 
perature and emissivity (Dunne et al. 2002), we calculate SFRs 
in excess of 1000 M yr -1 for redshifts past about 1. Of course 
the conversion between detected flux and inferred star formation 
rate is highly dependent on the dust SED, and can change by fac- 
tors of 10 for changes in temperature and /3 of only 2. Also, the 
simple relation between FIR luminosity and SFR may be differ- 
ent for these more luminous sources (Thr onson & Te lescolll986l 
iRowan-Robinson et all 19971) . 

Despite these uncertainties, it has been suggested (e.g. Blain 
et al. 2002, and references therein) that SCUBA sources can be as- 
sociated with the elliptical galaxies we see today via the following 
argument. Producing the local massive elliptical population with 
a homogeneous stellar distribution requires a sustained period of 
star formation on the order of 1000 M yr~ lasting about 1 Gyr. 
Based on results from Chapman et al. (2003) that place the bright 
SCUBA population at z > 1, the number of these galaxies per 
unit co-moving volume is comparable to the density of the local 
elliptical population. For example, if we take our estimate of the 
counts above 5 mjy and assume that they cover a redshift range be- 
tween 2 and 4 in a standard flat A-dominated model, we obtain a 
density of about 7 x 10 (ft Mpc) . These are thus rare and 
extremely luminous objects, with comparable number densities to 
galaxy clusters or quasars. 

If SCUBA sources really are associated with elliptical galax- 
ies, they should exhibit spatial clustering like their local counter- 
parts. There are other reasons one might expect detectable clus- 
tering; Extremely Red Objects (EROs) are very strongly clustered 
(Daddi etalJbood) . and seem to have a correlation with SCUBA 
sources. In general objects associated with major mergers should 
show high amplitude clustering (e.g. Percival et al. 2003). 



Although our analyses show no sign of clustering, the data are 
not powerful enough to rule it out. To improve on this we need more 
detected sources in order to bring down the Poisson error-bars. 
Also, the ERO and LBG clustering observations are taken from 
samples that exist at a common redshift (~ 1 in the case of EROs 
and ~ 3 for LBGs). Because of the strong negative K-correction, 
detected SCUBA sources are spread across a much wider redshift 
range, therefore diluting the clustering signal. Hence, progress can 
only be made with a larger survey (such as SHADES) that also has 
the ability to discriminate redshifts, even if only crudely. 

Although more studies are required to verify this claim, it is a 
reasonable hypothesis, and one with some testable predictions. Our 
new catalogue of SCUBA sources (Table 2) should allow for future 
detailed comparison with other wavelength data, which facilitate 
such tests. 



9 CONCLUSIONS 

This paper has presented the most complete accounting of sub-mm 
flux in the HDF-North region to date. We were able to demonstrates 
that > 4cr SCUBA detections are quite robust, being consistently 
detected in independent observations of the same area. Our cata- 
logue of sources was obtained using a careful statistical approach, 
involving simulations with the same noise properties as the real 
data. At 850 we were able to extract 19 sources above 4a and 
a further 15 likely sources above 3.5cr. Such a large list, in a field 
with so much multi-wavelength data, should be extremely useful 
for further studies. 

Our estimated source counts cover a wider flux density range 
than any other estimates, and given the careful completeness tests 
we carried out, they are likely to be more reliable than combining 
counts from different surveys. Our counts of SCUBA sources verify 
that significant evolution of the local LIRG population is required. 
Extrapolating a fit to these counts to below 2 mjy can reasonably 
recover the entire FIB at 850 fim. Clustering, although anticipated 
to be strong, was not detected in our map, due largely to the limited 
number of sources. Several hundred sub-mm sources with at least 
some redshift constraint will be required to detect the clustering 
unambiguously. 

The power of the SCUBA observations in the HDF-N lies not 
in the detection of objects per se, but rather for the ability to com- 
pare them with the plethora of existing and upcoming deep maps 
of this region at a wide variety of other wavelengths. Some of these 
comparisons will be the focus of paper II. 
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